clear;
exec('circshift.sci');
stacksize('max');

N = 512;

//%%%%%%%%%%%%%%%%% формирование сигналов %%%%%%%%%%%%%%%%%

 Xadc = csvRead('D:\Program Files\ADC_data\nois3_50dBm_GN22.csv');
//[x1, x2, x3, x4] = read_hex_data('D:\Program Files\ADC_data\nois3_60dBm_GN22\',1);
//Xadc = [x1; x2; x3; x4];

x = Xadc(:,1:$);//данные по строкам I-Q I-Q и т.д.
[num, len] = size(x);

//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

//%%%%%%%%%%%%%%%%%%%% MATRIX %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

main = [1 2];
E = x(main,:);
E_ = E(:,1:N);
//%  ???? ????????
//% Z = [x ; circshift(x',1)'; circshift(x',2)'; circshift(x',3)'; circshift(x',4)'; circshift(x',5)'];
//% Z(main,:) = [];
Z_data = x;
Z_data(main,:) = [];
n_delay = 4;
Z = Z_data;

for i = 1:n_delay
    Z = [Z; circshift(Z_data',i)'];
end

//Z_ = Z_data(:,1:N);
//
//for i = 1:2:n_delay
//    Z_ = [Z_ ; Z( (i*(num-2)+1) : (i+1)*(num-2), 1:N)];
//    if i == n_delay then
//        break;
//    end
//    tmp = circshift(Z_data(:,1:N)',(i+1))';
//    Z_ = [Z_ ; tmp];
//end

Z_ = Z(:,1:N);



Rxx = Z_*Z_';
Rxd = E_*Z_';
//pause
invRxx = inv(Rxx);
wm = Rxd*invRxx;
out = E - wm*Z;
db = 20*log10(stdev(E(1,:))/stdev(out(1,:)));
disp(db);
db = 20*log10(stdev(E(2,:))/stdev(out(2,:)));
disp(db);
disp(stdev(out(1,:)));
disp(stdev(out(2,:)));

//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

//%%%%%%%%%%% MATRIX + SUBBAND %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

num_bw = 4;
bp = 0:0.5/num_bw:0.5;
order = 2*num_bw+1;
wft = zeros(1,order);
wfm = zeros(1,256);
fr = zeros(1,256);
for i=1:num_bw
    [wft_t, wfm_t, fr_t]=wfir('bp',order,[bp(i) bp(i+1)], 're',[3 1.8]);
    wft(i,:) = wft_t;
    wfm(i,:) = wfm_t;
    fr(i,:) = fr_t;
end

y = zeros(num_bw*num, len);

for i=1:num_bw
    for j=0:num-1
        y(j*num_bw+i,:) = filter(wft(i,:), 1,x(j+1,:));//расположение данных по строкам 
//        %I1sub1 I1sub2 ... Q1sub1 Q1sub2 ... I2sub1 I2sub2 ... Q2sub1 Q2sub2 ...
    end
end

mainw = [ (main(1)-1)*num_bw+1:main(2)*num_bw ];

Ew = y(mainw,:);
Ew_ = Ew(:,1:N);

Zw_data = y;
Zw_data(mainw,:) = [];
n_delayw = 4;
Zw = Zw_data;
for i = 1:n_delayw
    Zw = [Zw ; circshift(Zw_data',-i)'];
end

Zw_ = Zw(:,1:N);
outw = zeros(2*num_bw,len);
for i=1:num_bw
    Zw__ = Zw_(i:num_bw:$,:);
//    for i = 1:n_delayw
//        Zw__ = [Zw__; circshift(Zw_(i:num_bw:$,:)')'];
//    end
    Ew__ = Ew_(i:num_bw:$,:);
    Rxxw = Zw__*Zw__';
//    invRxxw = inv(Rxxw);
    Rxdw = Zw__*Ew__';
    ww = Rxxw\Rxdw;
    outw(i:num_bw:$,:) = Ew(i:num_bw:$,:) - ww'*Zw(i:num_bw:$,:);
end

//Rxxw = Zw*Zw';
//invRxxw = inv(Rxxw);
//Rxdw = Zw*Ew';
//ww = invRxxw*Rxdw;
//outw = Ew - ww'*Zw;

dbw = 20*log10(stdev(sum(Ew(1:$/2,:),'r'))/stdev(sum(outw(1:$/2,:),'r')));
disp(dbw);
dbw = 20*log10(stdev(sum(Ew($/2+1:$,:),'r'))/stdev(sum(outw($/2+1:$,:),'r')));
disp(dbw);
disp(stdev(sum(outw(1:$/2,:),'r')));
disp(stdev(sum(outw($/2+1:$,:),'r')));
dbw_ = 20*log10(stdev(Ew(1:$/2,:),2)./stdev(outw(1:$/2,:),2));
// disp(dbw_);

// df = 1/(4*num_bw);
// f = df:2*df:0.5-df;
// bar(f,dbw_,1);
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



